import numpy as npimport pandas as pdimport osimport csvimport composite as ciimport inspectfrom functools importreduceimport plotly.express as px# from xml.etree import ElementTree as ET# import requests# import json# import webbrowser
Code
# for VSC users, Latex in Plotly is not working# Workaround from https://github.com/microsoft/vscode-jupyter/issues/8131import plotlyfrom IPython.display import display, HTMLplotly.offline.init_notebook_mode()display(HTML('<script type="text/javascript" async src="https://cdnjs.cloudflare.com/ajax/libs/mathjax/2.7.1/MathJax.js?config=TeX-MML-AM_SVG"></script>'))
Code
# show all output from cellsfrom IPython.core.interactiveshell import InteractiveShellInteractiveShell.ast_node_interactivity ="all"# last_exprpd.set_option('display.max_columns', 20)
Code
# create output folder if not thereifnot os.path.exists("../output"): os.makedirs("../output")
Country names dict
Code
# load countries dictionarywithopen("../data/countryAbbr.csv") as csv_file: reader = csv.reader(csv_file) country_to_abbrev =dict(reader)# invert the dictionaryabbrev_to_country =dict(map(reversed, country_to_abbrev.items()))# add additional abbreviationsabbrev_to_country["EL"] ="Greece"abbrev_to_country["GB"] ="United Kingdom"
Beach litter originating from national land-based sources that ends in the beach (%)
PERCENT
900
EN_MAR_BEALIT_BV
Beach litter originating from national land-based sources that ends in the beach (Tonnes)
TONNES
900
EN_MAR_BEALIT_EXP
Exported beach litter originating from national land-based sources (Tonnes)
TONNES
900
EN_MAR_BEALIT_OP
Beach litter originating from national land-based sources that ends in the ocean (%)
PERCENT
900
EN_MAR_BEALIT_OV
Beach litter originating from national land-based sources that ends in the ocean (Tonnes)
TONNES
900
EN_MAR_CHLANM
Chlorophyll-a anomaly, remote sensing (%)
PERCENT
2784
EN_MAR_CHLDEV
Chlorophyll-a deviations, remote sensing (%)
PERCENT
4131
EN_MAR_PLASDD
Floating plastic debris density (count per km2)
NUMBER
3
14.2.1
EN_SCP_ECSYBA
Number of countries using ecosystem-based approaches to managing marine areas (1 = YES; 0 = NO)
NUMBER
33
14.3.1
ER_OAW_MNACD
Average marine acidity (pH) measured at agreed suite of representative sampling stations
PH
903
14.4.1
ER_H2O_FWTL
Proportion of fish stocks within biologically sustainable levels (not overexploited) (%)
PERCENT
422
14.5.1
ER_MRN_MPA
Average proportion of Marine Key Biodiversity Areas (KBAs) covered by protected areas (%)
PERCENT
6049
14.6.1
ER_REG_UNFCIM
Progress by countries in the degree of implementation of international instruments aiming to combat illegal, unreported and unregulated fishing (level of implementation: 1 lowest to 5 highest)
NUMBER
882
14.7.1
EN_SCP_FSHGDP
Sustainable fisheries as a proportion of GDP
PERCENT
5880
14.a.1
ER_RDE_OSEX
National ocean science expenditure as a share of total research and development funding (%)
PERCENT
159
14.b.1
ER_REG_SSFRAR
Degree of application of a legal/regulatory/policy/institutional framework which recognizes and protects access rights for small-scale fisheries (level of implementation: 1 lowest to 5 highest)
NUMBER
882
14.c.1
ER_UNCLOS_IMPLE
Score for the implementation of UNCLOS and its two implementing agreements (%)
PERCENT
63
ER_UNCLOS_RATACC
Score for the ratification of and accession to UNCLOS and its two implementing agreements (%)
PERCENT
63
Check official data by indicator
Code
# check sdg indicator, change SeriesCode to the one you want to checksdg14check = sdg14[sdg14["GeoAreaName"].isin(countries)]sdg14check.loc[ sdg14check["GeoAreaName"].str.contains(r"^(?=.*United)(?=.*Kingdom)"), "GeoAreaName"] ="United Kingdom"sdg14check = sdg14check[sdg14check["SeriesCode"] =="EN_MAR_BEALITSQ"].pivot_table( columns="TimePeriod", index="GeoAreaName", values="Value", aggfunc="mean")sdg14checkmissingCountries(sdg14check)
TimePeriod
2015
2016
2017
2018
2019
2020
GeoAreaName
Belgium
1.520290e+07
2.704512e+05
2.273158e+05
1.915215e+04
6.645361e+05
NaN
Bulgaria
2.319362e+05
1.924616e+05
2.292929e+06
2.618936e+06
2.589610e+06
2.159091e+06
Croatia
4.538235e+07
4.292519e+05
NaN
NaN
5.534965e+04
NaN
Cyprus
3.804995e+05
1.230420e+06
4.879745e+06
4.901154e+06
NaN
1.879470e+07
Denmark
1.289859e+05
1.961387e+05
6.390451e+04
2.080706e+05
2.490390e+05
1.977920e+05
Estonia
NaN
NaN
1.058452e+06
NaN
NaN
NaN
Finland
NaN
NaN
2.172614e+04
1.652134e+05
NaN
NaN
France
6.693829e+06
4.550833e+06
1.906840e+05
1.066566e+06
1.329666e+06
6.046173e+05
Germany
1.010614e+05
2.106574e+06
4.334803e+04
3.861386e+05
1.695590e+04
NaN
Greece
2.866051e+06
1.806335e+06
1.039325e+06
1.432935e+06
1.515755e+04
7.201374e+04
Ireland
1.666727e+04
1.233589e+05
2.938770e+05
2.351199e+04
4.400202e+02
NaN
Italy
9.195292e+06
7.008116e+06
2.331030e+06
6.977531e+05
2.116006e+05
NaN
Malta
1.205128e+06
NaN
NaN
5.260832e+03
NaN
NaN
Netherlands
3.870849e+06
4.560344e+04
5.221688e+04
1.038472e+05
7.550317e+03
NaN
Poland
NaN
8.504996e+05
2.077893e+05
7.328124e+05
NaN
NaN
Portugal
3.096130e+06
1.237300e+06
2.573624e+05
1.537330e+06
1.605057e+05
1.029442e+04
Romania
1.460573e+05
2.776233e+06
1.341007e+07
7.194247e+05
NaN
NaN
Spain
1.531321e+06
4.088206e+06
2.193543e+06
1.436452e+06
1.257716e+06
2.463317e+05
Sweden
7.623596e+03
2.239108e+04
1.722629e+04
6.489057e+03
NaN
NaN
Missing countries:
Lithuania
Latvia
United Kingdom
14.1
(a) Index of coastal eutrophication
We use Gross Nitrogen Balance: Max-Min transformation where the max-value is the maximum value in the period 2012-2021 and the min-value is zero.
Proportion of national exclusive economic zones managed using ecosystem-based approaches
As of now, it is a binary variable: Number of countries using ecosystem-based approaches to manage marine areas
Data is scarce and the next round of data collection is in 2025. More information here
Spatial comparison within the EU is not expected to yield great differences and we neglect target 14.2.
14.3
Average marine acidity (pH) measured at agreed suite of representative sampling stations
We use two indicators:
Greenhouse gas emissions under the Effort Sharing Decision (ESD): Max-Min transformation where the max- and min-values are the maximum and minimum in the assessment period (2012-2021)
Carbon Emissions per capita: Max-Min Transformation where the max-value and min-values are the maximum and minimum in the period 2012-2021.
1. Greenhouse gas emissions under the Effort Sharing Decision (ESD)
Several sources (see code)
Code
# Greenhouse gas emissions under the Effort Sharing Decision (ESD), Million tonnes CO2 equivalent (Mt CO2 eq), European Environment Agency# https://www.eea.europa.eu/data-and-maps/data/esd-4ghgESD = pd.read_excel("../data/EEA_ESD-GHG_2022.xlsx", sheet_name=1, skiprows=11)ghgESD = ghgESD[ghgESD["Geographic entity"].isin(countries)]ghgESD = ghgESD.set_index("Geographic entity")ghgESD = ghgESD.dropna(axis=1, how="all")negativeValues(ghgESD)mr = ghgESD.columns[-1] # most recent yearghgESD = ghgESD *1_000_000
Code
# Build ESD allocation with two sources and interpolation for years previous to 2013# Allocation for 2020 target# https://ec.europa.eu/clima/ets/esdAllocations.do?languageCode=en&esdRegistry=-1&esdYear=&search=Search¤tSortSettings=allocation2020 = pd.read_xml("../data/esdAllocations2020.xml", xpath=".//ESDAllocations/ESDAllocationInfo")allocation2020["Country"] = ( allocation2020["ESDMemberState"] .map(abbrev_to_country) .fillna(allocation2020["ESDMemberState"]))allocation2020 = allocation2020[allocation2020["Country"].isin(countries)]allocation2020 = allocation2020.pivot_table( columns="ESDYear", index="Country", values="Allocation", aggfunc="mean")# Allocation for 2030 target# https://eur-lex.europa.eu/legal-content/EN/TXT/HTML/?uri=CELEX:32020D2126allocation2030 = pd.read_html("../data/esdAllocations2030.html")[13][1:]allocation2030.columns = allocation2030.iloc[0]allocation2030 = allocation2030[1:]allocation2030.loc[ allocation2030["Member State"].str.contains("Netherlands"), "Member State"] ="Netherlands"allocation2030 = allocation2030[allocation2030["Member State"].isin(countries)]allocation2030.set_index("Member State", inplace=True)for col in allocation2030.columns: allocation2030[col] = allocation2030[col].apply(lambda x: str(x).replace("\xa0", "") )allocation2030 = allocation2030.astype(int)# Merge 2005 values with 2020 and 2030 allocations. Interpolate for years 2006-2012allocationESD = ghgESD[[2005]].merge( allocation2020.merge( allocation2030, left_index=True, right_index=True, how="outer" ), left_index=True, right_index=True, how="outer",)allocationESD.columns = allocationESD.columns.astype(int)allocationESD[list(range(2006, 2013))] = np.nan # create empty columns for 2006-2012allocationESD = allocationESD[list(range(2005, 2031))]allocationESD.interpolate(axis=1, inplace=True)# Calculate score for ESDscoreESD =100- (100* (ghgESD - allocationESD) / allocationESD)scoreESD[scoreESD >100] =100scoreESD = scoreESD.ffill(axis=1) # fill empty values with last available yearscoreESD[[2013, 2016, 2021]]missingCountries(scoreESD)
2013
2016
2021
Geographic entity
Belgium
100.000000
99.633951
100.000000
Bulgaria
100.000000
100.000000
100.000000
Croatia
100.000000
100.000000
100.000000
Cyprus
100.000000
100.000000
86.945820
Denmark
100.000000
100.000000
100.000000
Estonia
100.000000
100.000000
100.000000
Finland
100.000000
96.549221
100.000000
France
100.000000
100.000000
100.000000
Germany
100.000000
99.619542
100.000000
Greece
100.000000
100.000000
100.000000
Ireland
100.000000
99.312005
96.146602
Italy
100.000000
100.000000
100.000000
Latvia
100.000000
100.000000
100.000000
Lithuania
100.000000
100.000000
100.000000
Malta
92.959862
85.673820
100.000000
Netherlands
100.000000
100.000000
100.000000
Poland
100.000000
99.344500
100.000000
Portugal
100.000000
100.000000
100.000000
Romania
100.000000
100.000000
100.000000
Spain
100.000000
100.000000
100.000000
Sweden
100.000000
100.000000
100.000000
United Kingdom
100.000000
100.000000
100.000000
No missing countries
Alternative (dif-ref) score for ESD (measuring this way does not make a lot of sense to me)
Code
# Alternative (dif-ref) score for ESD (measuring this way does not make a lot of sense to me)scoreESD2020 =100-100* ( ghgESD.loc[:, :2020].subtract(allocationESD[2020], axis=0)).divide(allocationESD[2020], axis=0)scoreESD2030 =100-100* ( ghgESD.loc[:, 2021:].subtract(allocationESD[2030], axis=0)).divide(allocationESD[2030], axis=0)scoreESD1 = scoreESD2020.merge(scoreESD2030, left_index=True, right_index=True)scoreESD1[scoreESD1 >100] =100# scoreESD1[[2012, 2018, 2021]]
2. Carbon emissions per capita
Several sources (see code)
Code
# Population on 1 January, Eurostat# https://ec.europa.eu/eurostat/databrowser/view/tps00001/pop = pd.read_csv("../data/tps00001.csv")pop["geo"] = pop["geo"].map(abbrev_to_country).fillna(pop["geo"])pop = pop[pop["geo"].isin(countries)]pop = pop.pivot_table( columns="TIME_PERIOD", index="geo", values="OBS_VALUE", aggfunc="mean")
Code
# CO2, TOTX4_MEMONIA, THS_T thousand tonnes, Eurostat# https://ec.europa.eu/eurostat/databrowser/view/ENV_AIR_GGE/co2 = pd.read_csv("../data/env_air_gge.csv")co2["geo"] = co2["geo"].map(abbrev_to_country).fillna(co2["geo"])co2 = co2[co2["geo"].isin(countries)]co2 = co2.pivot_table( columns="TIME_PERIOD", index="geo", values="OBS_VALUE", aggfunc="mean")mr = co2.columns[-1] # most recent yearco2pc = co2 / pop *1000*1000## tonnes co2 per capita# maxMin transformation as (max-x)/(max-min) * 100 --> converts to 0-100 scale with 0=worst, 100=bestco2pc = ( (co2pc.loc[:, 2012:mr].max().max() - co2pc.loc[:, 2012:mr])/ (co2pc.loc[:, 2012:mr].max().max() - co2pc.loc[:, 2012:mr].min().min())*100).round(2)co2pc[[2012, 2016, mr]]missingCountries(co2pc)
TIME_PERIOD
2012
2016
2021
geo
Belgium
42.05
45.74
53.10
Bulgaria
64.56
65.88
56.04
Croatia
78.36
80.32
79.02
Cyprus
25.37
21.33
28.96
Denmark
46.29
51.71
63.46
Estonia
17.51
9.15
41.36
Finland
54.35
52.19
57.25
France
71.19
74.19
78.78
Germany
44.53
47.30
58.87
Greece
46.28
58.48
59.09
Ireland
25.41
22.20
29.20
Italy
67.06
74.56
75.56
Latvia
86.16
79.85
69.73
Lithuania
80.85
79.24
74.57
Malta
47.73
79.86
81.87
Netherlands
30.39
32.16
47.20
Poland
47.12
47.37
44.51
Portugal
73.10
71.04
81.12
Romania
81.28
88.51
86.97
Spain
62.10
67.35
73.73
Sweden
92.50
94.04
98.95
United Kingdom
NaN
NaN
NaN
No missing countries
Code
# Data to corrobarate the ENV_AIR_GGE data. Very close numbers and we have for 2021# Fossil CO2 emissions by country (territorial), million tonnes of carbon per year (1MtC = 1 million tonne of carbon = 3.664 million tonnes of CO2)# https://globalcarbonbudget.org/carbonbudget/co2T = pd.read_excel("../data/National_Fossil_Carbon_Emissions_2022v1.0.xlsx", sheet_name=1, skiprows=11)co2T = co2T.Tco2T.columns = co2T.iloc[0, :]co2T.columns = co2T.columns.astype(int)co2T = co2T.rename_axis(index="geo", columns=None)co2T = co2T[co2T.index.isin(countries)]mr = co2T.columns[-1] # most recent yearco2T = co2T *3.664# convert from carbon to co2co2pc = co2T / pop *1000_000## tonnes co2 per capita# maxMin transformation as (max-x)/(max-min) * 100 --> converts to 0-100 scale with 0=worst, 100=bestco2pc = ( (co2pc.loc[:, 2012:mr].max().max() - co2pc.loc[:, 2012:mr])/ (co2pc.loc[:, 2012:mr].max().max() - co2pc.loc[:, 2012:mr].min().min())*100).round(2)co2pc = co2pc.ffill(axis=1) # fill empty values with last available yearco2pc[[2012, 2016, mr]]missingCountries(co2pc)
2012
2016
2021
geo
Belgium
47.69
51.35
55.72
Bulgaria
69.89
71.90
73.62
Croatia
87.80
89.00
88.46
Cyprus
54.53
52.28
54.06
Denmark
65.29
70.80
82.76
Estonia
12.53
13.24
59.32
Finland
45.77
53.09
68.23
France
79.37
83.14
87.31
Germany
40.23
43.51
57.15
Greece
56.07
69.72
81.01
Ireland
55.34
53.56
62.32
Italy
68.21
75.71
78.69
Latvia
94.41
94.54
93.06
Lithuania
85.17
86.47
83.59
Malta
70.46
100.00
97.23
Netherlands
42.63
43.76
57.51
Poland
53.21
53.55
52.37
Portugal
85.49
84.33
92.02
Romania
87.58
92.56
90.60
Spain
75.18
78.05
83.89
Sweden
83.74
88.28
96.28
United Kingdom
60.82
73.99
84.41
No missing countries
14.4
Proportion of fish stocks within biologically sustainable levels
We use two indicators:
FMSY/F: catch-weighted average
B/BMSY: catch-weighted average
1. FMSY/F
2. B/BMSY
14.5
Coverage of protected areas in relation to marine area
We consider two indicators:
Coverage of protected areas in relation to marine areas: Distance to Reference transformation where the max-value is set to 30% and the min-value is 0%
OHI ‘Biodiversity’ Index : No further transformation
1. Marine protected areas (% of territorial waters)
Degree of implementation of international instruments aiming to combat illegal, unreported and unregulated fishing
We use two indicators:
Fishing Subsidies relative to Landings: Max-Min Transformation where the max-value and min-value are the maximum and minimum in the assessment period (2012-2021)
Proportion of total research budget allocated to research in the field of marine technology
We use two indicators:
Official UNSD indicator ER_RDE_OSEX: Max-Min transformation where the max-value and min-value are the maximum and minimum in the assessment period (2009-2017) of all European countries (not restricted to the countries in the analysis).
SAD/TAC: Catch-weighted TAC relative to Scientific Advice
Note: ER_RDE_OSEX data comes from Global Ocean Science Report (GOSR) 2020, and goes from 2013 to 2017. Data for 2009-2012 data is available in the UNstats archive (download csv for 29-Mar-19)
# %%time# National ocean science expenditure as a share of total research and development funding (%), UNstats archive (years 2009-2012)# https://unstats.un.org/sdgs/indicators/database/archive# # read old data by chunks to speed loading and save 14.a.1 to separate file# iter_csv = pd.read_csv('./data/AllData_Onward_20190329.csv', iterator=True, chunksize=1000, encoding='cp1252')# oResearchOld = pd.concat([chunk[chunk.Indicator == '14.a.1'] for chunk in iter_csv])# oResearchOld.to_csv('./data/archive14a1.csv', index=False)oResearchOld = pd.read_csv("../data/archive14a1.csv")oResearchOld = oResearchOld.pivot( index="GeoAreaName", columns="TimePeriod", values="Value")
Code
# National ocean science expenditure as a share of total research and development funding (%), UNstats# https://unstats.un.org/sdgs/dataportal/database# read official data and merge with archive dataoResearch = sdg14[sdg14["SeriesCode"] =="ER_RDE_OSEX"].pivot_table( columns="TimePeriod", index="GeoAreaName", values="Value", aggfunc="mean")oResearch = oResearchOld.merge( oResearch, left_index=True, right_index=True, how="outer")# use all countries in EuropeoResearch = oResearch[oResearch.index.isin(allEurope)].dropna(how="all", axis=1)# fill nan of year 2013 from new report with old report and 2016 with 2017oResearch[2013] = oResearch["2013_y"].fillna(oResearch["2013_x"])# oResearch[2016] = oResearch[2016].fillna(oResearch[2017])oResearch = oResearch[list(range(2013, 2018))]# weighted by EEZ areaoResearch = oResearch.merge(eez, left_index=True, right_index=True, how="outer")for col in oResearch.drop("Area_km2", axis=1).columns: oResearch[col] = ( oResearch[col] * oResearch["Area_km2"] / (oResearch["Area_km2"].sum()) )# maxMin transformation as (x - min)/(max-min) * 100 --> converts to 0-100 scale with 0=worst, 100=bestoResearch = ( (oResearch.loc[:, 2013:2017] - oResearch.loc[:, 2013:2017].min().min())/ ( oResearch.loc[:, 2013:2017].max().max()- oResearch.loc[:, 2013:2017].min().min() )*100).round(2)# fill nan with mean of columnfor col in oResearch.columns: oResearch[col] = oResearch[col].fillna(oResearch[col].mean())oResearch = oResearch[oResearch.index.isin(countries)]oResearch[[2013, 2016, 2017]]missingCountries(oResearch)
2013
2016
2017
Belgium
0.020000
0.000000
0.000
Bulgaria
0.100000
0.050000
0.060
Croatia
7.850000
6.764444
12.894
Cyprus
3.858182
6.764444
12.894
Denmark
3.858182
6.764444
12.894
Estonia
3.858182
6.764444
12.894
Finland
1.210000
1.230000
1.210
France
13.300000
6.764444
8.860
Germany
0.840000
0.650000
0.630
Greece
3.858182
6.764444
12.894
Ireland
3.858182
6.764444
100.000
Italy
9.550000
7.950000
5.030
Latvia
3.858182
6.764444
12.894
Lithuania
3.858182
6.764444
12.894
Malta
3.858182
6.764444
12.894
Netherlands
0.740000
0.620000
0.690
Poland
0.210000
0.130000
0.120
Portugal
3.858182
36.750000
12.894
Romania
0.890000
6.764444
12.894
Spain
7.730000
13.500000
12.340
Sweden
3.858182
6.764444
12.894
United Kingdom
3.858182
6.764444
12.894
No missing countries
2. SAD/TAC
Code
# If a country fishes only on fish stocks where assignment of TAC follows scientific advice, it would score 100
14.b
Degree of application of a legal/regulatory/policy/institutional framework which recognizes and protects access rights for small‐scale fisheries
We use two indicators:
OHI Artisanal Fishing Opportunities Index: No further transformation
Percentage of Fish Species Threatened: No further transformation
# Percentage of Fish Species Threatened# Waiting to get IUCN API key# I would not use time series comparison because of this: https://www.iucnredlist.org/assessment/red-list-index# https://fishbase.mnhn.fr/country/CountryChecklist.php?c_code=056&vhabitat=threatened&csub_code=# Last updated 2019# We can extract the data and calculate the percentage# https://en.wikipedia.org/wiki/ISO_3166-1_numeric# https://data.worldbank.org/indicator/EN.FSH.THRD.NO# Data is only for 2018###################### https://stats.oecd.org/Index.aspx?DataSetCode=WILD_LIFE## Here we have by country updated most recent.# Missing countries: Bulgaria,Cyprus,Croatia,Ireland,Malta,Romania,United Kingdomthreatened = pd.read_csv("../data/WILD_LIFE_10062023112434707.csv")threatened = threatened[ (threatened.SPEC =="FISH_TOT") & (threatened.IUCN =="THREAT_PERCENT")]threatened = threatened[["Country", "Value"]].set_index("Country")threatened = threatened[threatened.index.isin(countries)].rename( columns={"Value": 2021})threatened[2021] =100- threatened[2021]threatenedmissingCountries(threatened)
2021
Country
Belgium
79.577
Greece
90.526
Netherlands
76.289
Poland
80.159
Portugal
66.154
Sweden
84.848
Estonia
91.262
Italy
92.668
Finland
86.364
Germany
75.127
Latvia
97.701
Lithuania
95.238
France
84.049
Spain
96.338
Denmark
96.639
Missing countries:
Bulgaria
Cyprus
Croatia
Ireland
Malta
Romania
United Kingdom
14.c
Number of countries making progress in ratifying, accepting and implementing through legal, policy and institutional frameworks, ocean-related instruments that implement international law, as reflected in the United Nations Convention on the Law of the Sea
We use two indicators:
Participation in agreements of the International Marine Organization (IMO Participation Rate):
Measures under the Marine Strategy Framework Directive:
1. Participation in agreements of the International Marine Organization
# dictionary for country codes used by the GISISwithopen("../data/gisisDict.csv") as csv_file: reader = csv.reader(csv_file) gisisDict =dict(reader)gisisDict = {v: k for k, v in gisisDict.items()}
Code
# Participation in agreements of the International Marine Organization# https://www.imo.org/en/About/Conventions/Pages/StatusOfConventions.aspx Excel file with current status of IMO conventions# We get the historical data from the GISIS database: https://gisis.imo.org/Public/ST/Ratification.aspx# You need to create account to access data.# I tried to scrape the data but I am getting errors with Selenium and bs4.# I downloaded the html manuallygisisCountries = {k: v for k, v in gisisDict.items() if k in countries}listIMO = []for v in gisisCountries.values(): link ="https://gisis.imo.org/Public/ST/Ratification.aspx?cid="+str(v) listIMO.append(link)# for i in listIMO:# webbrowser.open(i)
Code
imoRatedf = pd.DataFrame(columns=["Country", 2012, 2018, 2021])# loop thru the html files in the folder and extract the datafor i inrange(len(os.listdir("../data/treatiesIMO/"))): imoRate = pd.read_html("../data/treatiesIMO/Status of Treaties » Ratification of Treaties{}.html".format( i ) )[4]# get the country name from the html file name by checking if string is in the list of countriesfor country in countries:if country in imoRate["Treaty"][0]: countryIMO = country imoRate.columns = imoRate.iloc[1] imoRate = imoRate[2:]# new column with the year of accession and denounced imoRate["accession"] = ( imoRate["Date of entry into force in country"] .str.extract("^([^(]+)") .fillna("") ) imoRate["denounced"] = imoRate["Date of entry into force in country"].str.extract(".*\\:(.*)\\).*" ) imoRate[["accession", "denounced"]] = ( imoRate[["accession", "denounced"]] .apply(pd.DatetimeIndex) .apply(lambda x: x.dt.year) )# count the number of treaties each country accessioned and didn't denounced by 2012, 2018 and 2021for i in (2012, 2018, 2021): imoRate[str(i)] = np.where( (imoRate.accession < i)& ((imoRate.denounced > i) | (imoRate.denounced.isna())),1,0, ) imoCount = ( countryIMO, imoRate["2012"].sum(), imoRate["2018"].sum(), imoRate["2021"].sum(), ) imoRatedf.loc[len(imoRatedf), imoRatedf.columns] = imoCount# calculate total possible treaties, apply dif-ref and convert to percentagetotalIMO =len(imoRate.dropna(subset=["Date of entry into force in country"]))imoRatedf = imoRatedf.set_index("Country").sort_index()imoRatedf =100-100* (totalIMO - imoRatedf).divide(totalIMO)imoRatedf = imoRatedf.astype(np.float64)imoRatedf = imoRatedf.apply(pd.to_numeric)imoRatedfmissingCountries(imoRatedf)
2012
2018
2021
Country
Belgium
78.431373
84.313725
90.196078
Bulgaria
76.470588
80.392157
82.352941
Croatia
72.549020
74.509804
74.509804
Cyprus
68.627451
70.588235
72.549020
Denmark
80.392157
88.235294
92.156863
Estonia
76.470588
78.431373
82.352941
Finland
76.470588
86.274510
90.196078
France
84.313725
90.196078
96.078431
Germany
80.392157
88.235294
88.235294
Greece
80.392157
84.313725
84.313725
Ireland
66.666667
68.627451
68.627451
Italy
72.549020
74.509804
74.509804
Latvia
80.392157
80.392157
82.352941
Lithuania
58.823529
62.745098
64.705882
Malta
60.784314
66.666667
66.666667
Netherlands
84.313725
90.196078
92.156863
Poland
80.392157
84.313725
86.274510
Portugal
68.627451
78.431373
86.274510
Romania
58.823529
62.745098
64.705882
Spain
86.274510
88.235294
88.235294
Sweden
82.352941
92.156863
94.117647
United Kingdom
78.431373
78.431373
78.431373
No missing countries
Measures under the Marine Strategy Framework Directive
See code for notes
Code
# The barplot here: https://eur-lex.europa.eu/legal-content/EN/TXT/PDF/?uri=CELEX:52018DC0562&from=EN# Comes from https://eur-lex.europa.eu/legal-content/EN/TXT/PDF/?uri=CELEX:52018SC0393# The analogous report for 2020 is here https://environment.ec.europa.eu/system/files/2023-04/C_2023_2203_F1_COMMUNICATION_FROM_COMMISSION_EN_V5_P1_2532109.PDF# But the assessment is much shorter. They refer the reader to a JRC report:# https://publications.jrc.ec.europa.eu/repository/handle/JRC129363# That report assesses all the descriptors, but it cannot be compared to the previous assessment.# Moreover, the source code and data are not available.# Overall, it is hard to make an indicator for the measures taken against pressure indicators by the MS.# Countries report different measures and data is poor.
Indicators aggregation
Given our ratio-scale full comparable indicators,\(I_{it}\), meaningful aggregation of \(N\) indicators into a composite indicator \(CI_t\) is obtained according to social choice theory by applying a generalized mean:
\[CI_t(\alpha_{it},I_{it},\sigma) = \left(\sum^N_{i=1}\alpha_{it}I^{\frac{\sigma-1}{\sigma}}_{it}\right)^{\frac{\sigma}{\sigma-1}} \quad \text{for} \quad t = 2012, 2018, 2021 \text{(or most recent)}\]
with weights \(\alpha_{it} > 0\) and $0 ≤ ≤ $. The parameter \(\sigma\) is used to quantify the elasticity of substitution between the different indicators. High (low) values of \(\sigma\) imply good (poor) substitution possibilities which means that a high score in one indicator can (cannot) compensate a low score in another indicator. Consequently, high and low values of \(\sigma\) correspond to concepts of weak and strong sustainability, respectively. Depending on \(\sigma\), one can obtain a full class of specific function forms for the composite indicator.
We define:
\(\sigma_{Target} = 0.5\) and \(\sigma_{Target} = 10\)
Code
# fmt: off# missing 14.2, mortality, biomass, tacCatch, sadTac, msfdvarDf = [ nitro, wasteG, wasteR,scoreESD, co2pc, mpa, ohiBio, fseScore, ohiLive, ohiTour, oResearch, ohiArt, threatened, imoRatedf]varNames = [ "nitro", "wasteG", "wasteR", "scoreESD", "co2pc", "mpa", "ohiBio", "fseScore", "ohiLive", "ohiTour","oResearch", "ohiArt", "threatened", "imoRatedf",]# fmt: ondictIndicators =dict(zip(varNames, varDf))# stack variables in each dataframefor name, df in dictIndicators.items(): df = df.stack().to_frame().rename(columns={0: str(name)}) df.index.names = ["Country", "Year"] df.reset_index(inplace=True) df.Year = df.Year.astype(int) df.set_index(["Country", "Year"], inplace=True) dictIndicators[name] = df# merge all variables into one dataframe, forward and back fill by countryindicators = pd.concat(dictIndicators.values(), axis=1, join="outer")indicators = indicators.reset_index().sort_values(["Country", "Year"])indicators = indicators[indicators.Year.isin(list(range(2012, 2022)))]indicators = indicators.groupby(["Country"], group_keys=False).apply(lambda x: x.ffill().bfill())indicators = indicators.set_index(["Country", "Year"])indicatorsL = {"plastic": ["wasteG", "wasteR"],"14.1": ["nitro", "plastic"],# '14.2':[''],"14.3": ["scoreESD", "co2pc"],"14.4": ["mortality", "biomass"],"14.5": ["mpa", "ohiBio"],"14.6": ["fseScore", "tacCatch"],"14.7": ["ohiLive", "ohiTour"],"14.a": ["oResearch", "sadTac"],"14.b": ["ohiArt", "threatened"],"14.c": ["imoRatedf", "msfd"],}# calculate composite indicators for each target importing the function from composite.pytargets = indicators.copy()for target, indicator in indicatorsL.items():try: alpha =1/len(indicatorsL[target]) df = targets[indicator] sigma =10 targets[target] = ci.compositeDF(alpha, df, sigma)exceptKeyError:print("Missing indicator for", target)targets = targets[[i for i in targets if i.startswith("14")]]targets
Missing indicator for 14.4
Missing indicator for 14.6
Missing indicator for 14.a
Missing indicator for 14.c
%%time# calculate composite score for each target importing the function from composite.py# monte carlo for strong sustainability and one sigma for weak sustainabilityscoresStrong = ci.compositeMC(targets, years=[2012, 2016, 2021], simulations=10)scoresWeak = pd.DataFrame(ci.compositeDF(alpha =1/len(targets.columns), df = targets, sigma =10), columns=['scoreWeak'])
CPU times: user 228 ms, sys: 19.7 ms, total: 247 ms
Wall time: 242 ms
# merge all data and calculate min, max, and EEZ-weigthed average for indicators and targetsdata_frames = [indicators, targets, scoresStrong, scoresWeak]fullData =reduce(lambda left, right: pd.merge( left, right, left_index=True, right_index=True, how="inner" ), data_frames,)eezAverage = ( pd.DataFrame(fullData.stack().reset_index()) .merge(eez, left_on="Country", right_on="Territory", how="left") .rename(columns={0: "value", "level_2": "indicator"}))eezAverage["valueEEZ"] = eezAverage.value * eezAverage.Area_km2eezAverage = ( eezAverage.groupby(["Year", "indicator"]) .sum(numeric_only=True) .drop("value", axis=1))eezAverage["averageEEZ"] = eezAverage.valueEEZ / eezAverage.Area_km2
Code
# define function to plot boxplotsdef plotBox(df, xaxis_title=str, yaxis_title=str, xlegend=float, ylegend=float, maxShift=float, minShift=float): # start plots fig = px.box(df, x="indicator", y="value")# add Countries in the max and min# create annotation list, x is the x-axis position of the annotation annoList = [] maxCountriesD = {} x =0for s in df.indicator.unique():# get max value for indicator maxVal = np.max(df[df["indicator"] == s]["value"])# get countries with max value, if more than 4 countries, use * countries = df.loc[(df['value'] == maxVal)&(df['indicator']== s), 'Country'].valuesiflen(countries) >4: maxCountriesD[s] = countries countries ='*' countries =',<br>'.join(countries) annotation =dict( x=x,# y position is the max value y=maxVal, text= countries, yshift=maxShift, showarrow=False, ) annoList.append(annotation) x +=1 minCountriesD = {} x =0for s in df.indicator.unique():# get min value for indicator minVal = np.min(df[df["indicator"] == s]["value"])# get countries with min value, if more than 4 countries, use * and store in dictionary countries = df.loc[(df['value'] == minVal)&(df['indicator']== s), 'Country'].valuesiflen(countries) >4: minCountriesD[s] = countries countries ='*' countries =',<br>'.join(countries) annotation =dict( x=x,# y position is the min value y=minVal, text= countries, yshift=minShift, showarrow=False, ) annoList.append(annotation) x +=1# add EEZ-weigthed average values indicatorAverage = eezAverage.loc[(2021,df.indicator.unique()), :].reset_index()for s in indicatorAverage.indicator.unique():# get EEZ-weigthed average value for indicator averageVal =float(indicatorAverage[indicatorAverage["indicator"] == s]["averageEEZ"]) fig.add_scatter( x=[s],# y position is the average value y=[averageVal],type="scatter", mode ='markers', legendgroup ='EEZ average', name ='EEZ average', marker =dict(color ='black', size =6 )) fig.update_layout(annotations=annoList)# just to add the legend with one entry fig.update_traces(showlegend=False).add_scatter( x=[s], y=[averageVal],type="scatter", mode ='markers', legendgroup ='EEZ average', name ='EEZ average', marker =dict(color ='black', size =6 ))# change legend position, axis titles fig.update_layout( legend=dict( x=xlegend, y=ylegend, traceorder="normal", font=dict( family="sans-serif", size=12, color="black" ), bgcolor="White", ), xaxis_title=xaxis_title, yaxis_title=yaxis_title,# yaxis_range=[-20,100] ) fig.show()
# Get the targets for 2020 and 2030 in percentage# Member State greenhouse gas emission limits in 2020 and 2030 compared to 2005 greenhouse gas emissions levels# There targets for 2020 and for 2030# https://www.eea.europa.eu/data-and-maps/figures/national-progress-towards-greenhouse-gas# Official journals with the data can be found at (https://climate.ec.europa.eu/eu-action/effort-sharing-member-states-emission-targets_en)# limitESR = pd.read_excel('../data/targetsESR/FIG2-154307-CLIM058-v2-Data.xlsx', sheet_name=1, skiprows=19, header=1, skipfooter=32)# limitESR = limitESR.rename(columns={'Unnamed: 0':'Country', '(Resulting) ESR target 2030 (AR4)':'2030ESRtarget','ESR limit for 2020':'2020ESRtarget', '2005 ESD BJ':'2005Level'})# limitESR = limitESR[['Country', '2020ESRtarget','2030ESRtarget','2005Level']]# limitESR.set_index('Country', inplace=True)# limitESR = limitESR[limitESR.index.isin(countries)]# #UK is not in the dataset, we need to add from the official journal# for i in countries:# if i not in limitESR.index:# print(i, 'is not in the dataset')# limitESR
Code
# These data SHALL NOT BE USED. See reason on why ENV_AIR_GGE is preferable for the calculation:# (https://ec.europa.eu/eurostat/statistics-explained/index.php?title=Greenhouse_gas_emission_statistics_-_emission_inventories)# (https://ec.europa.eu/eurostat/statistics-explained/index.php?title=Greenhouse_gas_emission_statistics_-_air_emissions_accounts&oldid=551152#Greenhouse_gas_emissions)# CO2, KG_HAB, Eurostat# https://ec.europa.eu/eurostat/databrowser/view/ENV_AC_AINAH_R2/# co2 = pd.read_csv('../data/env_ac_ainah_r2.csv')# co2 = co2[co2['geo'].isin(countryAbb)]# co2['geo'] = co2['geo'].map(abbrev_to_country).fillna(co2['geo'])# co2 = co2.pivot_table(columns='TIME_PERIOD', index='geo', values='OBS_VALUE', aggfunc='mean')# mr = co2.columns[-1] # most recent year# # co2[[2012, 2016, mr]]/1000
Code
# From 14.6, using Eurostat data for landing values.# Even though the USD-EUR discrepancy does not affect the ratio we calculate,# we get today's exchange rate to convert landing values to USD and have a consistent unit# Solution source: https://stackoverflow.com/a/17250702/14534411# r = requests.get('http://www.ecb.int/stats/eurofxref/eurofxref-daily.xml', stream=True)# tree = ET.parse(r.raw)# root = tree.getroot()# namespaces = {'ex': 'http://www.ecb.int/vocabulary/2002-08-01/eurofxref'}# currency = 'USD'# match = root.find('.//ex:Cube[@currency="{}"]'.format(currency.upper()), namespaces=namespaces)# eurTOusd = float(match.attrib['rate'])# Landings of fishery products, TOTAL FISHERY PRODUCTS, Euro, Eurostat# https://ec.europa.eu/eurostat/databrowser/view/FISH_LD_MAIN/# landing = pd.read_csv('../data/fish_ld_main.csv')# landing['Country'] = landing.geo.map(abbrev_to_country).fillna(landing.geo)# landing = landing[landing.Country.isin(countries)]# landing = landing[['Country', 'TIME_PERIOD', 'OBS_VALUE', 'unit']]# landing['landingUSD'] = landing.OBS_VALUE * eurTOusd# landing.pivot_table(columns='TIME_PERIOD', index='Country', values='landingUSD', aggfunc='mean')